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ABSTRACT 

Observations indicate that outflows from massive young stars are more collimated during their early 
evolution compared to later stages. Our paper investigates various physical processes that impacts the 
outflow dynamics, i.e. its acceleration and collimation. We perform axisymmetric MHD simulations 
particularly considering the radiation pressure exerted by the star and the disk. We have modified the 
PLUTO code to include radiative forces in the line-driving approximation. We launch the outflow from 
the innermost disk region (r < 50 AU) by magneto-centrifugal acceleration. In order to disentangle 
MHD effects from radiative forces, we start the simulation in pure MHD, and later switch on the 
radiation force. We perform a parameter study considering different stellar masses (thus luminosity), 
magnetic flux, and line-force strength. For our reference simulation - assuming a 3OM star, we find 
substantial de-collimation of 35% due to radiation forces. The opening angle increases from 20° to 
32° for stellar masses from 2OM to 6OM . A small change in the line- force parameter a from 0.60 
to 0.55 changes the opening angle by ~ 8°. We find that it is mainly the stellar radiation which 
affects the jet dynamics. Unless the disk extends very close to the star, its is too small to have much 
impact. Essentially, our parameter runs with different stellar mass can be understood as a proxy for 
the time evolution of the star-outflow system. Thus, we have shown that when the stellar mass (thus 
luminosity) increases (with age), the outflows become less collimated. 

Subject headings: accretion, accretion disks - ISM: jets and outflows - methods: numerical - MHD - 
radiative transfer - stars: formation - stars: massive 



1. INTRODUCTION 

Outflows and jets are integral processes of star forma- 
tion. They are believed to be essential for the angular 
momentum evolution of the cloud core and the protostar 
- either directly as stellar winds or indirectly by chang- 
ing the structure and the evolution of the surrounding 
accretion disk. Also, outflows from young stars provide 
an important feedback mechanism to return mass and 
energy into the ambient medium from which the young 
star is born. 

Most of the current understanding about the forma- 
tion and propagation of jets/outflows comes from obser- 
vations of low-mass stars. In this case we are fortunate 
to know the leading dynamical parameters such as out- 
flow velocity, density, and temperature (e.g. Hartigan & 
Morse 2007, Ray et al. 2007). However, in case of mas- 
sive young stars many of these parameters are poorly 
known, although large mult i- wavelength studies for mas- 
sive star forming regions have been done over the last 
decade (Beuther et al. 2002a, Stanke et al. 2002, Zhang 
et al. 2005, Lopez-Sepulcre et al. 2009, Lopez-Sepulcre 
et al. 2010, Torrelles et al. 2011). These studies suggest 
that outflows are an ubiquitous phenomenon not only for 
low-mass stars, but also in massive star forming regions. 

The standard framework for the launching process of 
jets or outflows from low-mass protostars (and most 
probably also for extragalactic jets) is the model of a disk 
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wind accelerated and collimated by magneto-centrifugal 
and magnetohydrodynamic forces (Blandford & Payne 
1982, Pudritz et al. 2007). A number of numerical sim- 
ulations of jet formation have been performed which all 
confirm this picture of self-collimated MHD jets for low 
mass stars (Ouyed & Pudritz 1997, Krasnopolsky et al. 
1999, Fendt & Cemeljic 2002, Ouyed et al. 2003, Fendt 
2006, Fendt 2009). Its applicability has also been demon- 
strated for the case of extragalactic jets (e.g.Komissarov 
et al. 2007, Porth & Fendt 2010). 

The formation of the jets around young, still forming 
high-mass stars takes place in the deeply embedded cold 
dust and gas cocoons exhibiting large visual extinction 
of the order 100 to 1000 mag (Arce et al. 2007). Recent 
progress in observations of outflows from young high- 
mass stars comes from mm and cm wavelengths. Early 
low-spatial-resolution single-dish studies suggested that 
massive outflows may have a lower collimation degree 
than those of their low-mass counterparts (Shepherd & 
Churchwell 1996). Also, one of the famous outflow, the 
Orion-KL system, exhibits a more chaotic and not colli- 
mated structure (e.g., Schulz et al. 1995). However, later 
studies indicated that the collimation degree of jets from 
high-mass protostars can be as high as from low-mass 
regions (e.g. Beuther et al. 2002a, Beuther et al. 2002b, 
Gibb et al. 2003, Garay et al. 2003, Brooks et al. 2003, 
Beuther et al. 2004, Davis et al. 2004). Typical outflow 
rates in high-mass systems are estimated to be a few 
times 10 -5 to 10 -3 M yr _1 , implying accretion rates on 
the same order of magnitude (e.g., Beuther et al. 2002b). 

The originally puzzling result that different studies 
found different degrees of jet collimation for high-mass 
star-forming regions could qualitatively be resolved when 
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it was realized that these studies in fact targeted differ- 
ent evolutionary stages. While jets found in the youngest 
star forming regions appear similarly collimated as their 
low-mass counterparts, jets from more evolved (massive) 
stars exhibit much broader outflow cones. To account 
for this effect Beuther & Shepherd (2005) proposed an 
evolutionary picture for high-mass outflows, where the 
jet formation starts with similar MHD acceleration pro- 
cesses compared to the low-mass models during the ear- 
liest evolutionary stages. However, as soon as the central 
sources gain significant mass, other processes come into 
play, for example, the radiation from the central star or 
more turbulence at the base of the jet. Any combination 
of this and other processes was expected to lower the de- 
gree of collimation. However, the evolutionary sequence 
proposed by Beuther & Shepherd (2005) was largely ob- 
servationally motivated and could only provide a quali- 
tative explanation. A thorough theoretical treatment to 
account for this observational result is still missing - a 
topic which we now address in the present paper. 

How the radiation field affects the formation of a jet 
is not obvious a priori. In order to quantify and to 
disentangle the physical processes involved, a detailed 
numerical investigation is required. Essentially, stellar 
(and disk) radiative forces may affect jet acceleration and 
collimation directly (neglecting ionization, heating, and 
probably turbulent stirring for simplicity), but also in- 
directly by changing the physical conditions of the jet 
launching area, thus governing the mass loading or the 
initial entropy of the ejected jet material. For example, 
numerical MHD simulations have shown that jets with 
higher (turbulent) magnetic diffusivity are expected to 
be substantially less collimated (Fendt & Cemeljic 2002). 

Our previous studies (Vaidya et al. 2009) have shown 
that the inner accretion disk around massive protostars 
is ionized, has a high temperature and is gravitationally 
and thermally sufficiently stable in order to provide a 
suitable launching area for an outflow. Together with 
the recent observations of magnetic fields around high 
mass protostars (e.g Vlemmings 2008; Vlemmings et al. 
2010,Girart et al. 2009) , this indeed supports the picture 
of a scaled- up version of low-mass stellar jet formation. 

In this paper we will present a detailed investigation 
of how a strong radiation field impacts the structure and 
dynamics of a magneto- hydrodynamical driven jet. Mo- 
tivated by the presence of strong jets and outflows in 
massive star formation, we apply the standard picture of 
MHD jet formation known for astrophysical jets and put 
it in the physical environment of a massive young star. 

2. MODEL SETUP: JET FORMATION FROM MASSIVE 
YOUNG STARS 

In this section we discuss the model setup applied for 
our numerical study concerning the formation of jets and 
outflows from massive young stars. The central point is 
that we are going to consider the main features of the 
standard model of MHD jet formation which is well es- 
tablished for low-mass young stars or Active Galactic 
Nuclei also for high-mass young stars. Our model con- 
sists of the following essential ingredients (see Fig. 1). 

• A central massive young star which is rapidly evolv- 
ing in mass M*, luminosity L*, and radius R*. 

• A surrounding accretion disk with high accre- 




R* = 0.01 ...0.1 AU 

Figure 1. Sketch showing our model setup of inner regions around 
a massive young star. Several constituents are considered: The 
disk outflow (dashed arrows) is launched along the magnetic flux 
surfaces. The inner most launching point is denoted by -Rinjet- 
This outflow is subject to radiation forces (white arrow) from the 
star and (potentially) from an inner hot accretion disk. The stellar 
radius R* , could be as large as ~ 100 R© , while the stellar magnetic 
field structure and strength B* is rather uncertain. The location 
of an inner disk radius i?i n ,disk (if existent) is also not known. 

tion rate, estimated to be of the order of M ^ 
1(T 3 M yr- 1 . 

• A jet launching inner accretion disk. The extension 
of this area towards the star is not known and may 
depend on the existence of a strong stellar mag- 
netic field and stellar radiation pressure. Instead 
of introducing an inner disk radius i?i n ,disk, we will 
refer to an inner jet launching radius i^injet? which 
we presume to be between 0.1 and 1.0 AU, and 
to which the length scale of the simulation will be 
normalized Iq = it^njet- 

• A magnetic field around the protostellar object. 
Magnetic fields are essential for generating colli- 
mated high-speed outflows. Observational indica- 
tion for such fields around high-mass protostars ex- 
ists. 

• A strong radiation field of the high-mass young star 
which may influence accretion and ejection pro- 
cesses. In case when the accretion disk reaches 
down to radii close to the stellar surface (no gap 
as in case of low mass stars), also the high disk lu- 
minosity may play a role for the outflow dynamics. 

In the following we briefly discuss the observational 
and theoretical background of these constituents and fi- 
nally mention the limits of our model and possible model 
extensions (a more detailed discussion is provided before 
the summary). 

2.1. The central massive young star 
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It has been proposed that outflows from massive young 
stars may follow an evolutionary sequence such that out- 
flows tend to be more collimated and similar to jets from 
low-mass stars in the early stages of stellar evolution, 
whereas at later times the outflows are less collimated 
(Beuther & Shepherd 2005). Since this intrinsically cor- 
responds to an evolutionary sequence in stellar mass, we 
have investigated simulations with different central mass, 
ranging from 20 M to 60 M . 

A higher stellar mass automatically gives rise to a 
faster outflow (supposed the relative launching radius is 
the same), just because the outflow originates deeper in 
the gravitational well. 

2.2. The accretion disk and the jet launching area 

Jets from low mass stars are thought to be launched 
from less than 1 AU of the inner accretion disk 
(e.g. Anderson et al. 2003, Ray et al. 2007). In T Tauri 
stars, these regions could be well studied via NIR inter- 
ferometry (e.g. Akeson et al. 2000, 2005, Dullemond & 
Monnier 2010). However, to probe regions < 100 AU 
around a young high-mass star is difficult due to the 
large extinction. We had therefore studied this region 
via semi-analytic modeling in a previous paper (Vaidya 
et al. 2009) to get handle on the physical properties of 
such disks. We have shown that here the disk may reach 
high temperatures ~ 10 5 K leading to sublimation of most 
of the dust and ionizing the bulk of the material. 

The inner disk radius in case of low mass stars is usu- 
ally estimated assuming magnetic pressure balance with 
the accretion ram pressure. Typical values are ~ 3 — 5 R*. 
For young massive stars, such an estimate is not possible 
due to lack of knowledge on stellar magnetic fields during 
the formation stage. In addition, radiative force from the 
bright luminous massive star could also influence dynam- 
ics of the inner disk. Although our semi-analytic mod- 
eling indicated that the disk could extend right down to 
the central star, detailed 3D models with accurate radia- 
tive treatment are required to get more insight in these 
close- by regions (Vaidya et al. in prep). 

We have chosen the inner launching point at a dis- 
tance of Iq = 1 AU from the star (i.e. similar to low mass 
stars). A value of Iq < 0.1 AU would imply a high rota- 
tion speed and a deep potential well, resulting in a faster 
outflow. At the same time the jet launching part of the 
disk would be much hotter and possibly result in a higher 
contribution to the radiation force from these hot inner 
parts. On the other extreme, large jet launching radii 
lo > 10 AU will result in slow outflows, barely affected 
by stellar and disk radiation forces. 

The density po at the inner launching point is used for a 
physical scaling. We estimate po from the observed mass 
fluxes, which are typically of the order of 10 -3 — 1O _5 M 
yr _1 , providing po ~ 10 -13 — 10 -15 gem -3 (§ 5.3). 

2.3. The magnetic field 

For many years the role of magnetic fields in massive 
star formation was not really known. However, recent 
observations have detected relatively strong magnetic 
fields in massive star forming regions (Vlemmings 2008). 
Polarimetric observations of the hot massive molecular 
core HMC G31.41 have revealed a large-scale hourglass- 
shaped magnetic field configuration (Girart et al. 2009). 



Beuther et al. (2010) detected a magnetic field aligned 
with the molecular outflow via polarimetric CO emis- 
sion. Further observations have detected synchrotron 
emission from the proto-stellar jet HH 80-81, indicating 
a ~ 0.2 mG magnetic field in the jet knots while the stel- 
lar mass of ~ 10 M is in the range of massive stars 
(Carrasco-Gonzalez et al. 2010). 

Vlemmings et al. (2010) have measured magnetic field 
strengths using polarization by 6.7 GHz methanol masers 
around the massive protostar Cepheus A HW2 and de- 
rive a line of sight (l.o.s) magnetic field strength ~ 23 mG 
at a distance of 300-500 AU from the central star. The 
magnetic field strength is parametrized by the plasma 
beta, /?o, which is the ratio of the thermal gas pressure 
to the magnetic pressure at the inner launching point 
lo. Our simulations are so far limited to 1 < j3o < 10 
by numerical and physical reasons. This translates into 
field strengths at 1 AU of ~ 100 times weaker than esti- 
mated by conserving magnetic flux using observed values 
at 300 AU. 

Note that it is not only the field strength but also the 
field distribution which affects the jet formation process, 
as it was shown by (Fendt 2006) by MHD jet forma- 
tion considering a wide parameter set of magnetic field 
and mass flux distribution along the jet launching area. 
As a result, simulations applying a concentrated mag- 
netic flux profile tend to be less collimated. We apply 
as initial field distribution the standard potential field 
suggested by Ouyed & Pudritz (1997) (§ 4.2). A central 
stellar dipolar field (see Fig. 1) is not (yet) supported by 
observations. 

2.4. The stellar and disk luminosity 

The massive young star produces a substantial lumi- 
nosity, which is supposed to dynamically change the out- 
flow structure. The dependence of the radiation force on 
the stellar luminosity is parametrized by the Eddington 
ratio T e , defined as the ratio of the radiation force due 
to electron scattering to the central gravity (see Table 
1). The characteristic values of this dimensionless pa- 
rameter are obtained from the stellar evolution model 
of Hosokawa & Omukai (2009). In order to study the 
impact of radiation forces on the dynamics of outflows, 
we first launch a collimated jet from a disk wind via 
MHD forces. Then, after the pure MHD jet has achieved 
a steady state, we initiate the radiation forces, and then 
compare their impact on the MHD jet. Also, simulations 
in which the radiation forces are considered from the be- 
ginning (which are computationally much more expen- 
sive) end up with a final structure of the outflow similar 
to the one obtained from the step-by-step method. 

We prescribe the radiation force by the line-driving 
mechanism introduced by Castor, Abbott, & Klein 
(CAK). Such a force is parameterized by two physical 
parameters k and a. The value of k is proportional to 
the total number of lines. The quantity a can be con- 
sidered as a measure for the ratio of acceleration from 
optically thick lines to the total acceleration (Puis et al. 
2000). Depending on the selection of lines, for massive 
OB stars typical values obtained for k range from 0.4- 
0.6, while a ranges between 0.3 and 0.7 (Abbott 1982). 
The process of line driving has also been applied to cata- 
clysmic variables (CVs) (Feldmeier & Shlosman 1999), as 
well as to hot and luminous disks around Active Galac- 
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tic Nuclei (AGN) (Proga et al. 2000, Proga & Kallman 
2004). 

Proga (2003) carried out numerical simulations driv- 
ing winds from hot luminous magnetized accretion disks 
of AGNs, assuming spherical geometry, an isothermal 
equation of state, and an initially vertical magnetic field 
structure. Here we consider a potential field which is 
hour-glass shaped, and an adiabatic equation of state. 



2.5. Limitations of our model setup 

Our paper, for the first time, provides a quantitative 
study of the interplay between radiative and MHD forces 
on outflows launched from the vicinity of young massive 
stars, applying high-resolution axisymmetric numerical 
simulations. However, a few critical points can be raised 
which may limit the applicability of our model and which 
should probably be considered in forthcoming investiga- 
tions. 

From the general point of view there is the lack of true 
knowledge concerning a number of important parameters 
as discussed above. One important question is the loca- 
tion of the inner jet launching radius. If it is identical 
with an inner disk radius - where is that inner disk ra- 
dius located, if it exists at all? Is there a strong stellar 
magnetic field which could open up a gap between the 
stellar surface and the disk as it is known for low-mass 
young stars? 

Further, our disk model is taken as a boundary con- 
dition steady in time. As the star evolves, also the disk 
structure and accretion rate may evolve in time. This 
question could only be answered by simulations solving 
also for the disk structure. We do not explicitly incor- 
porate heating and cooling but simply consider the adi- 
abatic expansion and compression in the jet. 

Another question is the existence of a stellar wind. 
We know that OB stars have strong mass loss in form 
of stellar winds during later stages of their life times. 
The derived mass loss rates are typically of the order 
of 10 -6 M yr _1 . These winds are primarily radiation 
driven via the line-driving mechanism (e.g. Kudritzki 
& Puis 2000, Owocki 2009). The velocities derived are 
high, ranging from - 500 - 1000 km s" 1 and are usually 
supersonic. However, in case of high-mass young stars, 
no indication for such winds has been found so far possi- 
bly due to high obscuration. Nevertheless, a future study 
should implement the physical effect of a central stellar 
wind. 



3. BASIC EQUATIONS 

For our study, we carry out axisymmetric numerical 
ideal MHD simulations using the PLUTO code (Mignone 
et al. 2007). We have modified the original code to in- 
corporate source terms treating the line-driven forces 
from central star and disk, taking into account self- 
consistently the density and velocity distribution of the 
outflow. 

The MHD code considers the following set of equations. 
That is the conservation of the mass, momentum, and 
energy, 



dp 
dt 



+ (v- V)p + /A7-i7=0, 
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where p is the gas density, v the velocity vector, P the 
gas pressure, and B the magnetic field vector with the 
poloidal and toroidal components B p , B^. In order to in- 
clude the radiative forces F rad , the relevant source terms 
have been added in the momentum and energy conser- 
vation equation. The total energy density of the flow E 
comprises contributions from the internal energy e, the 
mechanical energy, and the magnetic energy, 

B (4) 
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The gas pressure in the flow is related to the density as- 
suming an adiabatic equation of state with the adiabatic 
index 7, 

P = (7 - l)pe. (5) 

The evolution of the magnetic field is governed by the 
induction equation, 



dB 

~dt 



= V x 



(V x if) . 



(6) 



We treat the ideal MHD equations without considering 
resistive terms. 

In addition to the above set of equations the code obeys 
the condition of divergence-free magnetic fields, V • B = 
0, using the constraint transport method. 

3.1. Prescription of radiation forces 

We do not explicitly consider radiative transfer, how- 
ever, we study the effects of momentum transfer by ra- 
diative forces on the outflow matter which is launched 
by MHD processes from the underlying disk. The total 
radiative force F rad comprises of four contributions - the 
acceleration due to continuum radiation from star / CO nt,* 
and disk / CO nt,disk> respectively, and, similarly, due to 

spectral lines from star /n ne ,* and disk /ii ne ,disk, respec- 
tively, 



TU'cld 



fa 



/< 



cont,disk 



line,* 



1 ine,disk« 



(7) 



In case of young massive stars, the stellar radiation is 
sub-Eddington. This immediately implies that the con- 
tinuum force do not contribute in modifying the dynam- 
ics of the MHD outflow as its strength is much below 
the gravitational pull of the central star. However, line 
forces could prove to be efficient in substantially enhanc- 
ing the continuum force by a so-called force multiplier, 
which is subject to complex theoretical studies of radia- 
tive transfer. This theory was developed first by Castor 
et al. (1975) who solved the radiative transfer equations 
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from spectral lines in moving atmospheres. According 
to their studies, the force due to line driving can be ex- 
pressed as a product of force due to continuum radiation 
and a force multiplier. 

The force multiplier depends on two parameters - k and 
a. A "general" parametrization for the force multiplier 
which is independent of ion thermal velocities v t h (see 
Eq A2 in the Appendix) has been introduced by Gayley 
(1995). In his formulation, the parameter /c, initially 
introduced by Castor et al. (1975), has been replaced by 
a parameter Q and the force multiplier can be expressed 
as follows, 



M(T) = kT~ a 



1 



\n • V(n • v) 
cr e cp 



(8) 



where T is the optical depth parameter that depends on 
the l.o.s velocity gradients (see eq. (A2)). The "strongest 
form" of this parameterization requires the ansatz Q = 
Qo> where Qo being the line strength of the strongest 
line. 

Typically for evolved massive stars, the parameter Q 
lies in the range of 1000-2000 (Gayley 1995). For our 
present study, we have applied the force multiplier con- 
sidering this "general" parameterization in its "strongest 
form". Thus, we have the two parameters Qo and a 
which define the force multiplier (Eq. 8) and hence the 
line forces. For simplicity, we assume constant values for 
these line force parameters for a particular simulation 
run. One of the important properties of the force mul- 
tiplier M(T) is that its values saturates to a maximum 
value M max when the medium becomes optically thin 
and the optical depth parameter approaches a minimum 
value, T < Tmin- The typical value for the maximum 
force multiplier is relatively constant, M max ~ 10 3 , de- 
pending only on the metallicity (Gayley 1995). 

The line force due to the central star is a product of 
the force due to continuum from the star and the force 
multiplier, 

flme,*=fcont,*M(T). (9) 

The most difficult part for the numerical realization of 
the line forces is to calculate the proper l.o.s velocity 
gradients \n-V(h-v)\ that appears in the formulation 
of force multiplier (Eq 8). Following the definition of 
this term by Rybicki & Hummer (1978), we express this 
gradient in terms of the rate-of-strain tensor e^-, 



n-V(n-v) = e^nji 



tr 2 v*j 



&0i 

dr\ 



ninj, (10) 



where in general v-^n, and ni are the components of ve- 
locity v, distance r, and the unit vector n, respectively. 
The different components of this tensor in cylindrical co- 
ordinates were calculated from Batchelor (1967). Proga 
et al. (1998) approximated the above sum to be equal to 
the most dominant term, corresponding to the radial gra- 
dient of flow velocity along the spherical radius. Simula- 
tions with a more accurate algorithm of including all the 
other terms to determine the l.o.s. velocity gradient us- 
ing the above equation have also been carried out (Proga 
et al. 1999). These simulations show that the qualitative 
features of winds are not changed as compared to the ap- 
proximate calculation of the gradients. Further, they are 



numerically very expensive. We can approximate that 
the region of consideration is far away from the star as 
the central object is a point source. Therefore, the gra- 
dient can simply set to be equivalent to cIVr/cIR, where 
R is the spherical radius and Vr is the gas velocity along 
the radius R. However, when using cylindrical coordi- 
nates instead of spherical coordinates we have to re- write 
equation (10) by transforming R to r, and thereby adding 
further terms, 



h • V(n • v) 



1 



dr 



dv r dv z 
dz dr 



>dv z 
dz 



(ii) 

where x = z/r. 

We have applied two of the line force components indi- 
vidually as source terms in order to disentangle their ef- 
fects on the dynamics of a (pure) MHD jet launched from 
the underlying disk. These two components are from the 
central star alone and other due to underlying hot ac- 
cretion disk. Disks around massive young stars accrete 
very rapidly with rates of 1O _5 M0 yr _1 to 1O _3 M yr _1 . 
Such high accretion rates imply very high accretion lumi- 
nosities, in particular in their very inner regions. In or- 
der to calculate the line forces due to the disk luminosity, 
proper geometric factors have to be taken into account. 
We apply a formulation similar to Pereyra et al. (2000) 
whose details are given in the appendix. 

4. NUMERICAL SETUP 
4.1. Grid setup & physical scaling 

We perform axisymmetric ideal MHD simulations of 
jet formation in the presence of radiative forces. Our 
simulations are carried out on a grid of physical size (r x 
z) = (52 1 x 152 Z ), where Iq denotes the physical length 
scale. The grid is divided into (512 x 1024) cells in radial 
and vertical direction, respectively. Within r < Iq and 
z < Iq the grid is uniform with a resolution of Sr = 5z = 
0.05, while for l < r < 50 / or l < z < 150 Zo tne 
grid is stretched with a ratio of 1.002739 and 1.001915 in 
radial or vertical direction, respectively. The remaining, 
very outer part of the grid is again uniform with cell size 
Sr — 0.125 and Sz = 0.20, respectively. 

Pure MHD simulations would be scale-free. However, 
when radiation is considered a physical scaling of the dy- 
namical variables becomes essential. Three scaling pa- 
rameters in physical units are used - the length scale Zq? 
the base flow density po and the Keplerian velocity vo, 
at the inner launching point (see § 2 for physical values 
used for scaling parameters.) 

All other quantities can be derived using the following 
definitions, 
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(12) 



The force multiplier defined in equation (8) can be ex- 
pressed in code units as follows 



M C (T C ) = 



1 - a \cr e cpolop c 



dv] 



dl 



(13) 



Quantities with subscript V are obtained from simu- 
lations in the dimensionless form where as the quantities 
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with subscript ; cgs ; are the ones required in the physi- 
cal units. We measure the time t c in units of Iq/vq. We 
denote the number of rotations at the inner launching 
point Iq as 7V rot = Iq/2ttvq = t c /2ii 

Using these normalizations the conservation of momen- 
tum (Eq. 2) can be written in dimensionless code units, 

+ (v c • \7 c )v c j = — ((V c x B c ) x B c ) 

Otc J PO 

-V c P c -p c V c ^ c + p c (F c rad ), (14) 

where the plasma beta, /3q = (8tt PoVq) / Bq , specifies the 
initial strength of the (poloidal) magnetic field at Iq and 
the radiation force 



T^rad 



Trad 



(15) 



Figure 2 displays the numerical setup used for our 
simulations. We also show the dynamically important 
forces that a fluid parcel experiences in an outflow from 
a young massive star - the gravitational force F grav it y 
by the central star, the Lorentz force components paral- 
lel ^Lorentz,|| (accelerating) and perpendicular F Lorentz? _L 
(collimating), the thermal pressure gradient F pressure , 

and the centrifugal force F cen t r if ug ai , which plays a vi- 
tal role in particular for accelerating the flow from the 
disk surface. The essential point of this paper is that 
we also consider radiation forces from the star and the 
disk (F ra diation, eq. (7)). We consider the most domi- 
nant radiative source terms which are those due to the 
line driving mechanism from the central massive star and 
the underlying (hot) disk (see Fig 1). 

4.2. Initial conditions 

We model the launching of the wind from an accretion 
disk representing the base of the outflow. The gravita- 
tional potential is that of point star located at a slight 
offset (— r g ,— z g ) from the origin to avoid singularity at 
r=z=0, 



$(r, z) oc ((r + r g ) 2 + (z + z g f 



-0.5 



(16) 



where r g = z g = 0.21 in all our simulations. 

Initially we prescribe a hydrostatic equilibrium with a 
density distribution 



p(r, z) oc ((r + r g ) 2 + (z + z g ) 2 ) 



-0.75 



(17) 



(Ouyed & Pudritz 1997). The thermal pressure follows 
an polytropic equation of state P = Kp 1 with 7 = 5/3. 
The magnitude of K = (7 — l)/7 is determined by the 
initial hydrostatic balance between gravity and pressure 
gradient. To begin with, all velocity components are van- 
ishing. Such an initial setup is typical for jet launching 
simulations (e.g. Ouyed & Pudritz 1997, Krasnopolsky 
et al. 1999, Fendt 2006). The initial hot hydrostatic dis- 
tribution of density and temperature does not play a sig- 
nificant role in determining the flow structure as it is 
eventually washed out by the cold [i.e. thermal pressure 
unimportant] MHD flow that is launched from the base. 

The initial magnetic field is purely poloidal (in the 
(r, z)-plane) with a distribution derived from the poten- 
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Figure 2. Schematic view of the numerical setup along with defi- 
nition of the boundaries conditions applied for most of the simula- 
tion runs. The lines indicate poloidal magnetic field lines anchored 
in the underlying accretion disk which is in slightly sub-Keplerian 
rotation. The dark black spot represents a fluid element frozen 
onto that field line. Vectors indicate and denote the various forces 
acting on the fluid element. See text for a detailed description of 
the different dynamically important forces. 

tial field B v = V x As with 



A, = V»- a + (* + *d) a -(*»+ £) (18) 

(Ouyed & Pudritz 1997), where za is considered as the 
dimensionless disk thickness. Such a magnetic field con- 
figuration is force free. Together with the hydrostatic 
density distribution this implies an initial setup in force 
equilibrium. The initial magnetic field strength is pre- 
scribed by choice of the plasma-beta /3q at the inner 
launching point. 

4.3. Boundary conditions 

To choose the correct physical boundary conditions is 
of utmost importance for numerical simulations as they 
describe the astrophysical system under consideration. 
In the present setup, we have to deal with four boundary 
regions (see Fig. 2). 

4.3.1. Axial boundary 

Along the jet axis an axisymmetric boundary condition 
is applied. The normal and the toroidal components of 
vector fields (i.e. B n: B^^v^vs) change sign across the 
boundary, whereas the axial components are continuous. 
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The density and the pressure are copied into the ghost 
zones from the domain. 

4.3.2. Equatorial Plane boundary 

The "disk boundary" is divided in two regions - the 
inner gap region extending from the axis to the inner 
launching radius, r < Zq? and further out the disk region, 
r > Zq? from where the outflow is launched. 

The setup of this boundary is the most crucial for the 
outflow simulation. This is an "inflow" boundary with 
the boundary values determining the inflow of gas and 
magnetic flux from the disk surface into the outflow. Spe- 
cial care has to be taken to consider the causal interaction 
between the gas flow in the domain and in the ghost cells. 

In order to have a causally consistent boundary condi- 
tion (Bogovalov 1997, Krasnopolsky et al. 1999, Porth 
& Fendt 2010), we impose the four physical quanti- 
ties p, ^ F ,P, E^. The toroidal electric field vanishes, 

= (v x B)^ — 0, as result of the ideal MHD ap- 
proximation. Thus, poloidal velocity and the poloidal 
magnetic field are parallel at the boundary, v p \\B p . The 
angular velocity of the field line (Ferraro's iso-rotation 
parameter) Q F , which is one of the conserved quanti- 
ties along the field line is fixed in time along the bound- 
ary. We have chosen a Keplerian profile along the disk 
Q (r) oc r -1,5 . The poloidal velocity at the boundary 
is "floating", i.e. copied from the first grid cell into the 
ghost cell each time step. Thus, the mass flux at the disk 
boundary is not fixed but consistently derived by causal 
interaction between the outflowing gas and the boundary 
value. 

The inner gap region is prescribed as hydrostatic pres- 
sure distribution. However, gas pressure gradient, grav- 
ity, and centrifugal acceleration is considered for the ra- 
dial force-balance in the disk, leading to a sub-Keplerian 
rotation 

^kep 

with x < 1- Solving for the radial balance delivers the 
density profile along the boundary 

Pdisk(r,z) = P^ z )' 

Therefore, there is a density contrast between disk and 
initial corona (the domain). In order to avoid numerical 
problems at the interface of between gap and "inner disk" 
(ie. the inner jet launching radius), we smoothen this 
transition using the Fermi function, considering 

= i X0 \ 

X Vl + exp(-10(r-l))y " 

The Fermi function is resolved by 16 grid cells. The 
pressure distribution along the boundary is fixed so as 
to have a cool disk which is three times denser than the 
initial corona. Thus, the wind emerging from such a 
cool disk has a similarly cool temperature (lower than 
the initial corona). Quantitatively, the average initial 
coronal temperature within the inner jet launching area 
is around 5 x 10 5 K, while the temperature of the disk 
surface (i.e., the boundary) in the same area reduces to 
5 x 10 4 K at 10AU. 



Table 1 

Summary of the dimensionless parameters to study the impact of 
radiation forces on the outflow dynamics 



Parameter Definition 



Eddington ratio ^ _ a e L* 

(proxy for stellar luminosity) e — 4ttcGM* 



Luminosity ratio = L acc _ GM*M 8 

(proxy for disk luminosity) ^ ~ L* i?*L* 



Radius ratio ^ _ R* 

(proxy for stellar radius) — lo 



Initial plasma beta „ _ SttPq 

(proxy for magnetic field strength) — Bq 



Line- force parameters Qo and a prescribes M(T) using Eq. (8) 



The magnetic field in the boundary is allowed to evolve 
in time to avoid any current sheets along the interface. 
When the simulation starts, the rotating disk winds-up 
the poloidal field and induces a toroidal field component. 
Constraining the field line angular velocity Q F to be con- 
stant in time, we need to adjust the rotational velocity 
of the gas in the boundary, 

^ = m F + (19) 

where the mass load of a field line r\ = (pv p )/B p is again 
a conserved quantity in stationary MHD. 

4.3.3. Outflow boundaries 

The right and top boundaries (see Fig. 2) are defined 
as outflow boundaries. The canonical outflow conditions 
(zero gradient across the boundary) are imposed to all 
scalar quantities and vector components, except for the 
magnetic field. The toroidal magnetic field component 
B^ in the boundary is obtained by requiring the gra- 
dient of the total electric current I to vanish, whereas 
the poloidal field components are estimated by setting 
the toroidal current density to zero (Porth & Fendt 
2010). Having this new outflow condition implemented 
in the standard PLUTO code, we ensure that there is no 
artificial collimation due to boundary effects. 

5. PARAMETER SURVEY 
5.1. Choice of governing parameters 

The three main parameters we apply for a compre- 
hensive study are the (i) central stellar mass M*, (ii) 
plasma-beta /3q, and (iii) the line force parameter a (see 
Table 1). We also find that the magnitude of the den- 
sity at the launching base of the flow po> which is a free 
parameter in our simulations, strongly affects the flow 
characteristics (see discussion in § 5.3). 

We have performed a large parameter study with re- 
spect to the central star. The stellar mass considered in 
our model ranges from 20 M to 60 M . The size and 
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Figure 3. Variation of the dimensionless parameter T e with stel- 
lar mass M* from Hosokawa &; Omukai (2009). The solid black line 
indicates the stellar mass evolution in the Hosokawa model. Stars 
represent the Eddington ratio T e for different stellar masses in our 
parameter survey § 6.1.5. The vertical dashed lines mark different 
evolutionary stages for massive young star accreting rapidly at the 
rateof lO-^Moyr -1 . The curve is obtained using tabulated evolu- 
tionary parameters kindly provided by Takashi Hosokawa (private 
communication) . 



the luminosity of these stars are obtained from the lit- 
erature stellar evolution models (Hosokawa & Omukai 
2009). The dimensionless parameter that controls the 
luminosity of the central star is T e (see Table 1). Its 
variation with the stellar masses considered in our model 
is shown in Fig. 3. (Also see § 6.1.5) 

We have carried out simulations for three different val- 
ues of magnetic field strength /3q = 1.0,3.0,5.0. We fix 
the density at the base of the flow to the same value as 
for the simulations for different stellar mass. Thus, for 
30 M star and an inner jet launching radius of 1 AU, the 
different /?o values correspond to a poloidal magnetic field 
ofll.5G, 6.6G, or 5.1G at the inner launching point, re- 
spectively. The detailed description of the different sim- 
ulations is shown in Table 2. The jet is launched by 
pure MHD processes. As it achieves a steady state af- 
ter say A/mhd inner rotations, we switch on the radiative 
line forces by the central star and run the simulation for 
another TVrad inner rotations. 

The stellar line force in our model is quantified by the 
parameters Qo and a. In principle, these parameters de- 
pend on the degree of ionization in the flow, however, for 
simplicity, we neglect this dependence and assume that 
their values remain the same for a particular simulation 
run. The line force is independent of Qo? however, the 
its dependence on a is strong Gayley (1995). We fix the 
line strength parameter Qo = 1400 for all runs. We carry 
out simulations with three values of a = 0.55,0.60,0.65. 
Although the differences in a seems to be small, they 
lead to considerable variation in the magnitude of the 
line forces (Gayley 1995). These values of a and Qq are 
consistent with the empirical values obtained for evolved 



massive stars (Abbott 1982, Gayley 1995). Since, we have 
no empirical values for these parameters during the for- 
mation stage, the present model uses similar values of 
line force parameters as obtained from evolved massive 
stars. 

5.2. Quantifying the degree of collimation 

There are in principal several options to quantify the 
collimation degree of jets. Fendt & Cemeljic (2002) sug- 
gested to compare the mass fluxes in radial and verti- 
cal direction of the jet as a measure of the jet collima- 
tion. However, in general the choice of a "floating" inflow 
boundary condition and the subsequent time-dependent 
change of mass fluxes makes it difficult to use them as 
a measure for the degree of collimation. Therefore, in 
this work, we quantify the jet collimation by the opening 
angle <p of the field lines (which are equivalent to velocity 
streamlines in a steady state). For comparison, we mea- 
sure the opening angle along a certain field line at two 
critical points along that field line viz. the Alfven point 
<pA and the fast magneto-sonic point 0p. Note that this 
is not the asymptotic jet opening angle. 

In order to estimate the amount of de-collimation 
A£[%] by radiation forces, we measure the angular sep- 
aration of field lines resulting from TVmhd disk rotations 
from those with after TVrad disk rotations, 



AC[%](s) = 100 



<t>{s, TVrad) ~ jV M HD; 

^mhd) 



(20) 



where s measures of path along the field line. A positive 
value of A£[%] implies de-collimation whereas negative 
values would imply collimation of the jet by radiative 
forces. 

5.3. The outflow density at the disk surface - po 

The radiative force term F c rad is a product of the force 
multiplier and the continuum force. The force multiplier 
(Eq. 8) depends on the physical scalings vq, lo and po 
(see Eq. 13). This does imply that each simulation run 
is unique to the chosen set of scaling parameters for a 
particular type of star. 

In case of massive stars these values are currently dif- 
ficult to measure very close to the star (see § 2). How- 
ever, observationally derived values of mass flow rates 
in molecular outflows and jets around massive stars 
(Beuther et al. 2002b, Zhang et al. 2005, Lopez-Sepulcre 
et al. 2009) can be used to constrain the densities with 
prior assumption of inner jet velocity. Typically mea- 
sured values of the mass outflow rate are of the order of 
10~ 3 to 1O~ 5 M yr" 1 . 

The density po at the base of the flow (z ~ 0) can be 
estimated by the mass flux launched from the base per 
unit time, 



M out = 27rpov ll 



rl /2 ~ q dr c 



(21) 



where the density at the base of the flow is p(r c , z c = 0) = 
por~ q . From Eq. (17), q = 3/2. We assume the matter 
is launched from the disk surface with Keplerian speed 

(i.e. v z (r c ,z c = 0) = vor c )• The physical scaling 
for the density, velocity and lengths are po 5 and Zq 
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Table 2 

Parameter study of stellar radiation line-force effects on MHD disk jets for different field strength and different stellar mass. We apply the 
same physical density at the jet base po = 5.0 x 10 -14 gcm -3 , the same inner launching radius Iq = 1 AU, and the same line-force 
parameter Qo = 1400 for all runs. The simulations are performed for A^mhd =319 inner rotations in pure MHD, followed by and 
iV^AD — 319 rotations, with switched-on radiative forces. A measure for the jet opening angle is given at the MHD critical points, </>a, 4>f 

and along the field line rooted at r = 5.0 AU. A£[%] denotes the percentage difference of the opening angles between the steady-state 
flows for pure MHD and including radiation force along the field line rooted at the same radius. The vertical mass flux M V ert[M0 yr _1 ] is 
measured along the top cells in the domain. The first ten simulation runs in the table apply a "floating" boundary condition for the 
injection velocity, while for last three runs a fixed-mass-flux boundary condition is applied (see § 6.2). 



Run ID 


M*[M ] 


(3o 


a 






max[AC[%]] 


Mvert 


[Moyr- 1 ] 


M30a055bl 


30 


1.0 


0.55 


21.98 


16.45 


5.0 


3.1 


X 


io- 


-5 


M30a055b3 


30 


3.0 


0.55 


25.63 


20.68 


17.0 


2.9 


X 


io- 


-5 


M30a055b5 


30 


5.0 


0.55 


26.05 


23.40 


34.5 


2.5 


X 


io- 


-5 


M30a055b5-c 


30 


5.0 


0.55 


25.55 


23.01 


33.9 


2.4 


X 


io- 


-5 


M20a055b5 


20 


5.0 


0.55 


20.44 


15.57 


5.8 


1.6 


X 


io- 


-5 


M25a055b5 


25 


5.0 


0.55 


23.27 


19.94 


21.9 


2.4 


X 


io- 


-5 


M45a055b5 


45 


5.0 


0.55 


26.84 


24.35 


37.7 


3.0 


X 


io- 


-5 


M50a055b5 


50 


5.0 


0.55 


28.86 


26.03 


42.0 


3.3 


X 


io- 


-5 


M60a055b5 


60 


5.0 


0.55 


31.45 


29.25 


56.0 


3.7 


X 


io- 


-5 


M60a060b5 


60 


5.0 


0.60 


23.73 


21.48 


28.7 


3.2 


X 


io- 


-5 


M60a065b5 


60 


5.0 


0.65 


20.93 


17.28 


11.7 


2.9 


X 


io- 


-5 


M25a055b3inj 


25 


3.0 


0.55 


45.84 


33.87 


30.2 


1.3 


X 


io- 


-5 


M30a055b3inj 


30 


3.0 


0.55 


47.71 


36.97 


37.2 


1.4 


X 


io- 


-5 


M50a055b3inj 


50 


3.0 


0.55 


48.76 


39.22 


42.1 


1.8 


X 


io- 


-5 



respectively, while r c is the non-dimensional length unit. 
Thus, the density scaling po can be expressed as 

-M-out i / l*c,max \ 

A> = o n ln • ( 22 ) 

27r ^0 L V r c,min/_ 

Using typical observed mass outflow rates for young mas- 
sive stellar jets, we calculate po ~ 10 -13 — 10 -15 gem -3 . 

The physical value of the density at the base of the flow 
is in the denominator of Eq. 13, and therefore affects the 
magnitude of the line-driving force significantly. Increas- 
ing the density by, say, a factor of 10, decreases the line 
force at least by a factor of 10". For the values of a 
considered her, this magnification could be as large as 3. 
Such a change in the radiative force has clearly a notable 
impact on the outflow dynamics, in particular on its de- 
gree of collimation. The description of our simulations 
with different base density is given in Table 3. 

6. RESULTS AND DISCUSSION 

The radiation field in massive star forming regions may 
play a crucial role in modifying the dynamics of outflows 
and jets. In order to disentangle the effects of radiative 
forces from the pure MHD jet formation, we decided to 
follow a two-step approach. We first i) launch a pure 
MHD disk jet and wait until it has reach a steady state 
(after about A/mhd = 320 rotations). We then ii) switch 
on the radiation line-forces and allow the jet flow to fur- 
ther develop into a new dynamic state. In some cases a 
new steady-state can be reached, in other cases an un- 
steady solution develops. 

For the radiation forces we consider line-forces exerted 
by i) the luminous massive young star, and those by ii) 
the surrounding accretion disk. 

6.1. Jet de- collimation by the stellar radiation field 



The radiation line-force from the central point star is 
determined by applying Eq. (9). For a detailed physical 
analysis of the effects of line- forces from the star alone, 
we have chosen simulation M30a055b5 as reference sim- 
ulation run (Table 2). The reference run is parametrized 
for a stellar mass of 3OM , surrounded by an accretion 
disk with an inner jet launching radius of 1 AU, and a 
density of 5.0 x 10 -14 gcm -3 at this radius. The initial 
poloidal magnetic field strength is fixed by a plasma- 
Po — 5.0. The radiative forces from the star are defined 
by the line force parameters Qo = 1400.0 and a = 0.55. 

The spatial distribution of the force multiplier M(T) 
as well as the specific stellar line force for the reference 
simulation are both shown in Fig. 4. The magnitude of 
the force multiplier peaks in the top-right low-density re- 
gions of the flow. This is expected as the force multiplier 
increases with decreasing density as shown in Eq.9. In 
order to calculate the true radiation force, the force mul- 
tiplier must be convolved with the continuum radiation 
from the star, resulting in a line-force distribution with 
a maximum close to the central star (see Fig. 4). 

The time evolution of the emerging jet^ is shown in 
Fig. 5, where we display the vertical jet velocity v z (r, z) 
and the poloidal magnetic field structure for reference 
simulation M30a055b5. We clearly see the change from 
the pure MHD flow (top panels) to the situation when 
stellar radiation forces are considered (bottom panels). 

Material which is injected from the disk surface (bottom 
boundary) is "frozen" on these field lines (ideal MHD). 
Initially, the plasma is accelerated magneto-centrifugally 
and gains substantial speed, producing a bow shock as 
it propagates. The bow shock leaves the domain after 
~ 60 rotations. Inertial forces of the outflowing mass 
flux induce a strong toroidal magnetic field component 

t Full movies are available for download at http://www.mpia. 
de/homes/f endt/vaidya/ research_bv . html 
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Table 3 

Parameter study of stellar radiation line-force effects on MHD disk jets for different jet density at the base of the flow with a stellar mass 
of 30 Mq. For all runs, plasma-/3o = 5, Iq = 1 AU and line-force parameters Qq = 1400 and a = 0.55. 



Runs p [gcm- 3 ] cf) A (p F max[AC[%]] M ver t [M© yr" 1 ] 

M30a055b5dl 3.0 x lO" 14 30.30 27.79 47.7 1.5 x lO" 5 

M30a055b5 5.0 x lO" 14 26.02 23.36 34.5 2.4 x lO" 5 

M30a055b5d2 1.0 x lO" 13 22.91 19.47 19.6 4.3 x lO" 5 

M30a055b5d3 5.0 x lO" 13 20.42 15.03 4.0 1.9 x 10~ 4 



Table 4 

Parameter study of disk radiation line-force effects on MHD disk jets. All simulation runs apply the same stellar parameters and 
line- force parameters as M30a055b5, however a different inner launching radius of only 0.1 AU. Additionally, two more dimensionless 
parameters for the disk radiation force fi = 0.4644 and A = 0.4969 are prescribed (see Table 1). 



Run 


p [gcm 3 ] 


max[AC[%]] 


M V ert [Moyr- 1 ] 


Comments 


Diskl 


5.0 x 10- 14 


-1.8 


6.1 x 10- 7 


reaches a steady state 


Disk2 


5.0 x 10~ 15 


3.1 


7.9 x 10~ 8 


remains unsteady 




10 20 30 40 50 10 20 30 40 50 

r fAUl r [AU1 



Figure 4. Contours of force multiplier M(T) and the specific line 
radiation force of a central star [in physical units] for the reference 
simulation M30a055b5 after Af ro t = 638 inner rotations. 

resulting in magnetic hoop stresses which self-collimate 
the magnetic field structure together with the hydrody- 
namic mass flux. Eventually, when all the dynamical 
forces along the field line are balanced again, the flow 
achieves a steady- state. The steady state magnetic field 
configuration obtained after the pure MHD flow is shown 
as white dashed lines in Fig. 5. 

At this stage we " switch on" the line-forces of the cen- 
tral star. Immediately, the emergence of a fast axial flow 
is visible. This flow is unsteady, first forming a knotty 
structure, and then, when approaching the upper bound- 
ary, stabilizes to a steady axial flow. Jet de-collimation 
by stellar radiation forces is indicated by the fact that 
the red solid field lines in the bottom panels of Fig. 5 do 
open-up significantly as compared to field lines in steady- 



state pure MHD simulation. 

When the radiative force is switched on, a shock front 
begins to propagate into the steady MHD jet (Fig. 5). 
This happens, because the additional radiation forces 
lead to an initial disturbance at the base of the flow, 
which is then propagated outwards. The effect of line- 
forces resulting in a series of propagating shocks is best 
seen in the poloidal velocity profile along a magnetic field 
line (see Fig. 6). The series of shocks eventually prop- 
agate out of the domain, and the flow attains a steady 
state again. However, the asymptotic outflow velocity 
which is then achieved is enhanced by a factor 1.5-2 
compared to the pure MHD flow. 

We have also carried out a run, M30a055b5-c, that in- 
cludes the radiation force right from the beginning along 
with the MHD flow. This flow evolves into the same 
configuration as our reference run where radiation force 
is "switched on" later. This proves that the initial con- 
ditions do not affect the final state of the flow. In our 
study, we prefer to use the step-by-step approach as it is 
computationally less expensive. 

In summary, we find that radiation forces modify the 
MHD disk jet essentially in terms of collimation, but also 
in terms of acceleration and terminal speed. 

6.1.1. Analysis of the force-balance in the outflow 

Here we investigate the main forces affecting the out- 
flow dynamics, comparing the magnitude of various force 
terms calculated at each grid point along the field line 
rooted at r = 5.0 AU for simulation run M60a055b5. Fig- 
ure 7 shows such a comparison of specific forces pro- 
jected parallel and perpendicular to the field line before 
and after considering line-forces due to stellar radiation. 

With respect to acceleration, the most striking feature 
is the enhanced specific pressure gradient force. When 
radiative forces are considered, the steady MHD flow ma- 
terial at higher altitudes is adiabatically compressed from 
the underlying accelerated material leading into a shock 
that propagates out of the domain. The resultant flow 
has higher temperature or an increase in thermal pres- 
sure. We also find that the gradient of the thermal pres- 
sure higher above in the flow may increase substantially, 
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Figure 5. Evolutionary sequence of the reference simulation M30a055b5. Shown is the jet vertical velocity distribution (color coding) 
and the poloidal magnetic field lines. The color bar represents the velocity scale in kms -1 for each row of panels. The series of images 
are taken at subsequent inner rotations as mentioned on the top of each image. The solid red lines are the poloidal magnetic field lines. 
To illustrate the impact of radiation forces, we also show the field lines obtained after the steady-state MHD flow for comparison as white 
dashed lines in the lower four panels. 



such that thermal pressure force and Lorentz force may 
become comparable. In terms of collimation, it is evident 
from Fig. 7 that the profile of the perpendicular compo- 
nent of the pressure gradient force along the field line be- 
comes flatter than in the pure MHD flow case. Thus, the 
pressure force which was not important for pure MHD 
flows, now becomes a significant factor in governing the 
dynamics of the initial flow acceleration (disk wind). 

Contrary, when the line-driving force is switched-on 
in the MHD jet, the centrifugal force becomes reduced 



by an order of magnitude particularly at high altitudes 
in the flow. In order to comprehend this variation of 
centrifugal force, it is useful to compare the conserved 
quantities of MHD, in particular the angular velocity of 
the field line Sip- In the pure MHD simulation is 
conserved throughout the simulation. This remains true 
also with radiative forces included, simply because 
is fixed as boundary value (see § 4.3). Therefore, the 
azimuthal velocity at the base of the flow is given by 
Eq. (19). At the same time, the additional line forces 
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Figure 6. Poloidal velocity along the field line rooted at r = 3 
AU for the two simulations M30a055b5 top panel and Disk2 bottom 
panel. The distance along the field line s is measured in AU. The 
various colored solid curves indicate the velocity profile at different 
inner disk rotations N TO t (blue for N TO t= 325, green for iV ro t = 333, 
red for 7V r ot=340, cyan for N YO t= 510, and magenta for N TO t = 
638). The black dotted line in both panels correspond to the last 
time step for the pure MHD flow. The solid yellow line shown in 
the bottom panel indicates the poloidal velocity profile after 1500 
inner rotations. 

accelerate the flow (close to the base) to higher poloidal 
velocities. The magnetic flux at the base is fixed as a 
boundary condition. Hence, the ratio of mass load to 
density r]/p in Eq. (19) increases, and since the toroidal 
magnetic field is negative, also the azimuthal veloc- 
ity decreases. This eventually leads to a decrease of the 
specific centrifugal force when line-forces are considered. 
The decrease in centrifugal force further implies that the 
outflow rotates slower, and that the acceleration close to 
the base is not only controlled by magneto-centrifugal 
forces, but also by the thermal pressure force and the 
radiation force (see Fig. 7). 

Close to the base of the flow, Lorentz forces are dy- 
namically not important. However, they peak at the 
Alfven point of that particular field line. Beyond the 
Alfven point, the Lorentz force becomes important and 
competes with other forces to govern the dynamics of the 
flow. When radiative force are considered, the Lorentz 
force close to the base of the flow and also at higher alti- 
tudes do not show significant deviations from its values 
in the pure MHD case. 

Line-driven radiative forces have considerable impact 
on the poloidal electric current distribution and also 
on the critical MHD surfaces in the jet. The contour 
plots shown in Fig. 8 compare the total poloidal current 
/ = r£>0 = J jpdA and position of critical surfaces be- 
fore and after considering line-driven forces in the refer- 
ence simulation M30a055b5. The figure shows that by 
adding radiative forces (stellar luminosity), the corre- 
sponding poloidal electric current density contours are 
shifted closer to the base of the flow, implying a lower 
toroidal field strength in these jets. We understand 
this results as a consequence of lack of jet rotation, thus 



less effective induction of toroidal magnetic field. 

Jet collimation in the conventional Blandford-Payne 
picture is caused by magnetic hoop stresses (-B^/r) of 
the azimutal magnetic field. Whereas, the larger toroidal 
magnetic pressure gradient aids in de-collimating the 
outflow. A lower would not only weaken the hoop 
stress but also reduce the magnetic pressure (~ BV). 
The decrease of toroidal magnetic pressure would lead to 
a lower magnetic pressure gradient force. So, a balance 
arises between the hoop stress that collimates the flow 
and the toroidal pressure gradient that de-collimates it. 
We observe that the field lines de-collimate to a wider 
configuration on addition of radiative forces, implying 
that the decrease of toroidal magnetic field has more im- 
pact on reducing of hoop stresses than reducing the mag- 
netic pressure gradient, thus eventually de-collimating 
the flow. 

Considering radiative forces also results in lower- 
ing the location of MHD critical surfaces (viz. slow 
magneto-sonic, fast magneto-sonic and Alfven surface) 
(see Fig. 8). The pure MHD flow is launched marginally 
super-slow and sub-Alfvenic. The flow near the axis is, 
however, sub-slow, due to the boundary condition of con- 
serving the initial hydrostatic density and pressure dis- 
tribution in the gap region between axis and disk. The 
flow speed at the critical points depends on the magnetic 
flux and mass density in the flow. These quantities re- 
main approximately unaltered when radiative forces are 
considered, thus the magneto-sonic wave speed remains 
similar as well. Since now the radiative line-forces ac- 
celerate the wind further to higher velocities, the critical 
surfaces shift to lower altitudes in the flow (i.e. close to 
the base of the flow). 

We finally note that radiative forces have both a direct 
as well as an indirect effect in modifying the collimation 
properties of the flow. The radiation force from the star 
directly affects the flow dynamics at its base and close to 
the star simply by transfer of momentum (Eq. 2), even- 
tually leading to a flow de-collimation. Indirectly, the 
radiation field also enhances the thermal pressure in the 
flow (Eq. 3) such that the specific pressure force becomes 
comparable to the Lorentz force - further de-collimating 
the flow. The resulting opening angles as measured from 
our reference run at various critical points along the field 
lines are listed in Table 2. For M30a055b5, the fi- 
nal state of the jet has a steady mass outflow rate of 
2.5 x 10" 5 M o yr- 1 . The opening angles are 26° and 23° 
at Alfven and fast point respectively for the field line 
which has a foot point r = 5 AU. The maximum per- 
centage separation A£[%] between the opening angles of 
this field line in pure MHD flow and that for the flow with 
radiative forces is 34.5%. A significantly positive value 
of A£[%] indicates flow de-collimation. In the following 
sections, we describe the effects of various physical pa- 
rameters that could affect the dynamics of the outflow. 

6.1.2. Radiation field and magnetic field strength 

Here we discuss the interrelation between radiative 
forces and a variation of magnetic flux and their com- 
bined effect on outflow collimation. We compare sim- 
ulation runs M30a055bl, M30a055b3, and M30a055b5, 
all assuming the same stellar mass 3OM , an inner jet 
launching radius of 1 AU, and a density at the base of 
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Figure 7. Comparison of all specific force terms projected parallel[\eft] and perpendicular[right] to the field line rooted at r=5.0AU for 
the run M60a055b5. Colors specify different specific force terms in dyne cm 3 g — 1 - the Lorentz force (green), the gas pressure gradient 
(red), and the centrifugal force (blue). The dashed lines corresponds to the respective force terms for the pure MHD run M60a055b5 (at 
steady-state), while the solid lines represents the forces for the final time step of the simulation including radiative forces. The brown solid 
line represents stellar radiation force term for the final time step. The black and magenta vertical solid lines mark the position of the 
Alfven and the fast magneto-sonic surface for the steady state flow including radiative forces, whereas the corresponding dashed lines are 
for the pure MHD flow. 



the flow of 5.0 x 10 14 gcm 3 , while the initial magnetic 
field strength parametrized by /?o ranges from 5.1 G to 
11.5 G. 

The (initial) pure MHD runs exhibit characteristic dif- 
ferences. A lower f3o provides faster jets which are more 
collimated, simply because the larger field strength pro- 
vides larger Lorentz forces to collimate and also accel- 
erate the flow more efficiently. A lower field strength 
implies a lower magneto-sonic wave speed, resulting in 
magneto-sonic surfaces located closer to the disk surface. 
Thus, we cannot simply compare the opening angles at 
the critical surfaces as their positions vary in the pure 
MHD runs with different field strengths. Instead, we 
quantify the actual change of opening angle with and 
without considering radiative force A£[%] (i.e. in per- 
centage), as profile along the field line. 

This measure is shown in Fig. 9 for a field line rooted 
at r = 5AU. We see that A£[%] is significantly posi- 
tive, indicating that radiative forces do considerably de- 
collimate the flow. However, the effect of de-collimation 
is enhanced when the magnetic field strength is reduced 
(i.e. an increase of /3q). The asymptotic value of A£[%] 
along a field line changes from 35% to 5% when the initial 
magnetic field strength is increased by a factor of two. 

The above analysis for different magnetic field strength 
suggests that in the strong field case (/3q = 1) jet colli- 



mation is controlled by the magnetic forces alone. How- 
ever, when we decrease the field strength (from 11.5 G to 
5.1 G as normalized at 1 AU) the stellar radiation line- 
force begins to compete with the magnetic forces. The 
resulting jet has a much wider field structure and outflow 
opening angle as it is de-collimated as compared to its 
pure MHD counterpart (see Fig. 5). We conclude that 
radiative forces from the central star (for the chosen pa- 
rameters) will dominate the magnetic effects in control- 
ling the flow collimation for field strengths < 5 G (at 
1 AU). In addition to the field strength the field profile is 
also important in dynamical evolution of the jet flow (see 
§ 2.3). Observational estimates of the field strength in 
these close by regions of massive young stars will further 
help to narrow the results. 

6.1.3. Impact of the line-force parameter a 

The stellar line-force is a non-linear function of a. Gay- 
ley (1995) showed that a significant change in the mass 
flux can be obtained with very slight change in a, indicat- 
ing that this parameter is very sensitive for calculation 
of radiative line- forces. 

Simulation runs M60a055b5, M60a060b5, and 
M60a065b5 apply the same magnetic field strength 
(/?o = 5.0), but differ in the line- force parameter from 
a = 0.55 — 0.65. A first result is that the increase of 
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Figure 8. Contours of the total poloidal electric current distribu- 
tion (left), and the location of the MHD critical surfaces (right). 
In both panels, the pure MHD flow is represented by black solid 
lines, while the radiative MHD flow is shown with red solid lines. 
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Figure 9. Flow collimation along the field line rooted at r=5 
AU. Shown is the profile of the percentage separation of field lines 
between the pure MHD flow and the radiative MHD flow [Af[%]] 
as measure of flow collimation. The solid, dashed and dot dashed 
lines are for three different plasma-/3o = 5.0,3.0, 1.0, respectively. 

a from 0.55 to 0.60 decreases the mass outflow rate 
significantly by ~ 16%. We will discuss this interrelation 
in § 6.2. However, also the outflow opening angles are 
affected as indicated by the top panel of Fig. 10. The 
outflows with lower a become more de-collimated than 
the others with higher a. 

The above trend indicates that the efficiency of the line 



Figure 10. Jet collimation, jet density and radiation force. 
Shown is the opening angle (i.e. the field inclination) at the Alfven 
point (stars) and the fast point (dots) of the field line rooted at 
r=5 AU for runs with different a (top panel), and different jet base 
density po (bottom panel). The dotted and solid lines corresponds 
to the opening angle at the Alfven point, and the fast point of the 
same field line after the pure MHD flow with an initial plasma-/3o 
= 5.0. 

driving mechanism is increased with lower a. Physically 
a lower a implies a higher contribution of optically thin 
lines in accelerating the flow. Thus, the dominance of less 
saturated (i.e. less self-shadowed) lines in accelerating 
the flow results in an efficient line driving. 

In summary, we find that the runs with lower a have 
higher mass flux and are less collimated. Our results 
from studying different a are in qualitative agreement 
with the analytical results from Gayley (1995). 

6.1.4. Impact of the density po 

In this section, we describe how collimation and ac- 
celeration in the outflow are altered with the density at 
the base of the flow. As mentioned in § 5.3, the physical 
mass density applied to scale the numerical simulation is 
a free parameter in our setup, were we estimate its value 
from the observations using Eq. 22. 

Simulation run M30a055b5dl has the lowest density 
among all four runs, po = 3 x 10 -14 gcm -3 For this run, 
the opening angles measured at critical MHD surfaces are 
larger compared to our reference run. Quantitatively, the 
opening angles at the Alfven and the fast surface are 30° 
and 28°, respectively. Also, the maximum percentage 
change in field line opening angle A£[%] increases from 
35% for the reference run to a staggering high value of 
48% for field line rooted at r = 5 AU. In case of a high 
jet base density po, the change in field line opening angle 
A£[%] reduces to 4%, with the opening angles at the 
Alfven and fast surface being 20° and 15°, respectively 
(Fig. 10). 

These results correspond to an inverse correlation of 
the density at the base of the flow with the force multi- 
plier M(T) (see Eq. 13). Higher densities result in an op- 
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tically thick environment, thus increasing the the optical 
depth parameter T, and reducing the force multiplier. 
Thus, the radiative line- driving force approaches the 
limit corresponding to the continuum radiation force. Es- 
sentially, the continuum force for typical massive young 
stars is weak compared to other dynamically important 
forces (for e.g. gravity). 

In summary, we observe the outflow density as one of 
the leading parameters to govern the dynamics of the 
outflow by radiation forces. However, observationally 
this density is very difficult, even impossible to deter- 
mine. We thus have to rely on certain assumptions, 
or may apply estimates of the outflow flux to calculate 
this density. Our parameter survey considers density val- 
ues which are consistent with the observed mass fluxes. 
From our studies, we observe in general that for densities 
po < 10 -13 gcm -3 the line-driving force from the central 
star is very efficient in accelerating and de-collimating 
the outflow. However, a denser environment will dilute 
the influence of the line-driving mechanism. 

6.1.5. Stellar mass evolution and outflow collimation 

Motivated by the hypothesis of Beuther & Shepherd 
(2005) we have studied the outflow dynamics and col- 
limation for a sequence of different stellar masses (viz. 
from 20 M to 60 M ). The change in stellar mass im- 
plies a change in the central luminosity. Fig. 11 shows 
the variation of the opening angle at the MHD critical 
points in simulation runs with central stellar mass. It 
can be seen that the opening angle varies from 20° to 
32° for stellar masses from 20 M to 60 M . The curve 
rises linearly up 30 M , but then the opening angle does 
not change considerably for increasing mass. 

For the physical parameters of the stellar evolution 
such as stellar luminosity and radius, we have applied 
the bloating star model of Hosokawa & Omukai (2009) 
assuming an accretion rate of 1.0 x lO _3 M yr _1 (see 
fig. 3). In this model for massive young stars, the stel- 
lar luminosity increases rapidly from 6M to 30 M . A 
star with high accretion rate undergoes swelling and is 
considered to bloat up to 100 R . It is the entropy dis- 
tribution in these stars which causes them to expand 
in size that much. For stars < 10 M , a thin outer 
layer absorbs most of the entropy from its deep interiors, 
thereby increasing the entropy in this layers. The rapid 
increase of entropy in the outer layer of the star causes 
the star to increase in size (Hosokawa & Omukai 2009). 
When the stellar mass reaches ~ 10 M , the star then 
begins to shrink (Kelvin-Helmholtz contraction) until its 
mass reaches ~ 30 M . After that, the star approaches 
the main sequence and then follows the typical main se- 
quence evolution. The stellar luminosity does not change 
considerably from M* = 30 M to 40 M - in fact there is 
a slight decrease. Contrary, for stars with mass > 40 M 
the stellar luminosity increases with mass. This is re- 
flected as a slight dip in the variation of T e with mass 
between 30 M and 45 M . Thus, the values of T e de- 
rived for these two masses do not show considerable dif- 
ference (0.2381 and 0.2269 respectively). We conclude 
that the radiation force and, hence, also the dynamics 
of the outflow do not change to a great extent for these 
masses. 

For simplicity we keep the inner jet launching radius 
and the density at the base of the flow fixed for all mass 




30 40 50 

Stellar Mass : M* [M©] 



Figure 11. Same as figure 10, but for simulation runs for different 
stellar mass. 

runs. The radiation force is not only altered by changes 
in T e , but also due to the physical scaling that appears 
in the prescription of the force multiplier. In particular, 
this is the Keplerian velocity at the inner launching ra- 
dius Zq? which is natu rally different for different masses 
and scales with \/M*. One would then expect that the 
radiative force changes by a factor equivalent to \/M*. 
The effect of such a change with increasing central mass 
not only enhances the mass flux launched from the disk 
boundary but also imparts the resulting outflow a wider 
morphology (see Table 2). 

In summary, with the above parameter study of differ- 
ent central stellar masses, we observe that as the lumi- 
nosity (or mass) of the central star increases, the stellar 
radiation force becomes relatively stronger. This leads to 
a higher degree of outflow de-collimation confirming the 
above mentioned observational picture of outflow evolu- 
tion in massive young stars. 

6.2. Radiation field and jet mass flux 

The simulations presented so far have been performed 
with "floating" boundary conditions at the base of the 
flow (see § 4.3). While guaranteeing a causally correct 
boundary condition ab initio, the disadvantage is that 
the mass flux emerging from the underlying disk bound- 
ary cannot be prescribed a priori. In fact, the mass flux 
is self-consistently calculated each time step by ensuring 
continuity of outgoing waves between the domain and 
the ghost cells. Thus, in our approach applied so far we 
fix the initial flow density and float the vertical outflow 
velocity at the boundary. 

The total mass outflow rate (M ra d + M vert ) in physical 
units is shown in Fig. 12 for runs with different stellar 
mass and applying a flo atin g boundary condition. The 
total mass flux has a \/M* dependence for pure MHD 
runs on converting it to physical units. While for a 
steady-state radiative MHD flow, the physical mass flux 
also depends on the Eddington parameter T e , which is 
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Figure 12. Outflow mass flux for simulations applying a different 
stellar mass. The top panel shows the total mass outflow rates for 
each star in case of a pure MHD flow (stars), and for the same out- 
flows including radiative forces (dots). The percentage difference 
between these mass fluxes is shown in the bottom panel. 

related to the stellar luminosity. 

The percentage change between these mass outflow 
rates is shown in the bottom panel of Fig. 12. Inter- 
estingly, the profile of this curves is similar to the varia- 
tion of the opening angles (see Fig. 11), indicating that 
the increased mass flux could in fact play a role for de- 
collimation. Thus, the combination of both the "float- 
ing" boundary conditions and the disturbance of the gas 
physics at the base of the flow by radiative source terms 
leads to the enhanced mass flux, which modifies the col- 
limation of the flow. 

However, also the direct effect of radiative force on the 
flow can physically deflect the flow and eventually lead to 
de-collimation. Thus, the problem becomes quite com- 
plex as the direct influence of radiative force on the flow 
is mingled with its second order effects for e.g. increase 
in the mass flux. In order to single out the influence of 
direct de-collimation effects by radiative forces, we have 
performed simulations where we have fixed the outflow 
mass flux as a boundary condition. 

It is essential to inject the outflow with a velocity 
which is supersonic to begin with. We fix the mass flux 
at the boundary by setting the inflow velocity to v z = 
0.24 x ^Kep? resulting in a flow with slow magneto-sonic 
Mach number close to unity and indicating a slightly su- 
personic flow. In Fig. 13 we show the mass outflow flow 
rates derived along the top (z = z m ) and right (r = r m ) 
boundaries (in red and blue, respectively). We have cal- 
culated the radial and vertical mass outflow rates using 
following expressions, 




M vert = / rpv z \ Zm dr. (24) 

JO 

In the simulations with floating boundary conditions, 
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Figure 13. The variation of the outflow mass flux with time in 
simulation runs M30a055b3inj (top panel) and M30a055b5 (bottom 
panel). The red and blue lines show the variation of the vertical 
and radial mass fluxes, respectively (see Eqs. 23,24). The black line 
corresponds to the mass flux injected from the underlying disk. The 
green line is the total mass flux leaving the computational domain. 
Coinciding green and black lines proves conservation of the mass 
flux. 

we observe that the mass flow rate reaches a steady-state, 
first derived self-consistently in the pure MHD flow, and 
then suddenly increased on adding the radiative source 
term. Essentially, this rise is not seen when we apply 
fixed-mass-flux boundary conditions. However, the fact 
that the curves of radial and vertical mass flux intersect 
after ~ 60 disk rotations after switching-on the radiation 
force, is a clear signature of de-collimation. Thus, the 
modification in terms of collimation due to the direct 
impact of radiative forces on a flow injected with constant 
mass flux can be now understood with our simulation 
runs with a priori fixed- mass-flux boundary conditions. 

We have thus performed three simulation runs with dif- 
ferent stellar mass applying the above-mentioned mass- 
injection boundary condition (see Table 2). The most 
evident measure of de-collimation that can be observed 
is the ratio of vertical to radial mass flow rate. A higher 
ratio indicates a higher degree of jet collimation. This 
mass flux ratio for the chosen stellar masses of 25, 30, 
and 50 M changes from 1.12 after the steady- state pure 
MHD flow to 0.92, 0.85, and 0.80 respectively, when ra- 
diative force is considered. Note that the radiative forces 
also accelerate the flow, increasing the maximum veloc- 
ity for M30a055b3inj in steady-state radiative MHD by 
a factor of ~ 1.5 compared to pure MHD flow. 

We also find, that de-collimation is mainly due to the 
direct influence of the stellar line-driven force on the out- 
flow, as the amount of mass flux launched into the do- 
main is fixed for these runs (it is thus not a 2nd-order 
effect resulting from a higher mass flux due triggered by 
the radiation force). 

In summary, we conclude here that the line-driving 
radiation force from the star plays a significant role in 
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directly modifying the dynamics of previously launched 
MHD flow. This force not only aids in acceleration of 
the flow, but also pushes the outflow material away from 
the axis resulting in a substantially wider opening angle. 

6.3. Acceleration by the disk radiation field 

We have also carried out a number of simulations that 
involve only radiation forces caused by a (hot) accretion 
disk. In principle, in addition to the intrinsic accretion 
luminosity, the disk radiative force can be enhanced by 
considering irradiation from the central star (Proga et al. 
1998, Drew et al. 1998). For the present model, however, 
we do not consider irradiation. 

Here we present results of two of our simulations, de- 
noted as Diskl and Disk2 (see Table 4). We have ap- 
plied the parameter set as for the reference simulation 
M30a055b5, however, two more parameters governing 
the dimensionless accretion rate that appears in the ac- 
cretion luminosity (Eq. B5) have to be prescribed. That 
is the ratio \i of disk to stellar luminosity, and the ratio A 
of stellar radius to inner launching radius (see Table 1). 
Naturally, a disk extending to radii closer to the center 
of gravity would have much higher temperatures thus lu- 
minosity (see Vaidya et al. 2009 for an application to 
massive young stars) and the resulting radiation forces 
would be able to affect the outflow more. 

For the number values applied for these additional 
parameters we again follow the bloating star model 
(Hosokawa & Omukai 2009) assuming an accretion rate 
of 10 -3 M©yr -1 , and an inner jet launching radius of 
l = 0.1 AU, thus, obtaining fi = 0.4644 and A = 0.4969. 
Note that the smaller inner jet launching radius also im- 
plies a inner disk radius closer to the star, compared to 
our models discussed above. 

For simulation run Disk2, we have increased the magni- 
tude of the disk radiative force by decreasing the density 
at the inner jet launching radius by a factor of 10, thus 
assuming po = 5 x 10 -15 gcm -3 . As the radiation force is 
very sensitive to the density, the lower density at the jet 
base does increase the radiation force from the disk by 
a factor of three. Note that, physically, this increase of 
the disk radiation force could mimic the effects of stellar 
irradiation. 

We quantify the change in collimation degree in runs 
Diskl and Disk2 again with the parameter A£[%] (see 
Eq. 20). Simulation run Diskl has a maximum of 
AC[%] = -1.8, while run Disk2 has AC[%] = 3.1%. 
These number values are rather low compared to our 
simulations with a stellar radiation force indicating that 
the disk radiation force alone is not strong enough to de- 
collimate the flow. 

Interestingly, simulation run Disk2 shows some un- 
steady behavior close to the axis - a feature which is 
absent when only stellar radiation forces were consid- 
ered. This is consistent with the findings of Proga et al. 
(1999), who detect a similarly unsteady behavior in their 
simulations for the case when the radiative force is dom- 
inated by the underlying disk. The fact that we do not 
see this feature in simulation Diskl suggests that the ra- 
diative force in run Diskl is not comparable to the other 
forces that control the flow dynamics, and that this out- 
flow does more correspond to a pure MHD flow. The 
unsteady flow behavior in run Disk2 is also reflected in 
the poloidal velocity evolution as steadily propagating 



'wiggles' in the velocity profile along the field line (or 
streamline) v p (s) (Figure 6). 

As maximum outflow velocity for the simulations Diskl 
and Disk2 we obtain relatively high velocities of ~ 350 
kms -1 . This is mainly due to the above mentioned fact 
that the outflow is launched deeper in the potential well 
and as close as Iq — 0.1 AU from the central star. Com- 
pared to the pure MHD run, we see, however, that the jet 
which is affected by disk radiation forces (only) achieves 
slightly higher asymptotic velocities, as the MHD reaches 
only ^ 300 km s -1 . Thus, the disk radiation force affects 
the outflow primarily by slightly accelerating it. We do 
not see indication that the outflow collimation is affected, 
which is understandable since the disk radiation force 
acts mainly in vertical direction. 

In summary, jet acceleration and collimation is rather 
weakly affected by the disk radiative forces as their mag- 
nitude is orders of magnitude smaller compared to ra- 
diation forces from the central star. However, for other 
astrophysical jet sources such as CVs and AGNs, the ra- 
diative force from the disk could play a significant role. 
Applying the scaling relations and comparing the amount 
of energy radiated per unit area from the standard thin 
disk from a typical young massive star (MS) and a stan- 
dard white dwarf (WD) we obtain from equation Bl, 

D MS = n nfi ( 10- 3 M Q yT^ \ / 3OM \ / 0.1AIA " 3 
d wd U - UD ^ ioi6g S -i J \ 1M J Vl0 9 cm7 

(25) 

. This indicates that the efficiency of line driving due to 
disk force alone for such massive proto-stellar disks is less 
than the compact object accretion disks like CVs which 
have been the focus of previous studies by Pereyra & 
Kallman (2003). This ratio is even smaller if we consider 
disks around AGNs as considered by Proga & Kallman 
(2004). 

6.4. Limitations of our model approach 

The main goal of our study was to disentangle the ef- 
fects of radiative forces from the young star-disk system 
on a pure MHD outflow launched within the standard 
picture of magneto-centrifugal, magnetohydrodynamic 
jet formation. In this section we discuss the limitations 
of our model approach to the subject of jet formation 
from massive young stars. 

6.4.1. Prescription of radiation force 

The prescription of radiation force used for the present 
model does not explicitly consider the radiation transfer. 
Instead it implements the radiative force due to lines as a 
source term in the momentum and the energy equations 
(see § 3). This greatly simplifies the numerical task and 
modeling becomes computationally inexpensive. 

6.4.2. Possible existence of a stellar wind 

Our model does not consider a stellar wind from the 
massive young star. Observed outflows around young 
massive stars are usually thought to have velocities of 
about 200 - 500 kms -1 , a value much smaller than the 
speed of stellar winds in the evolved stages of OB stars (of 
about ~ 1000 kms -1 ). This difference in velocity scale 
could give a clue to difference in environment around the 
star. During the formation stage, the star is surrounded 
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by dense gas and dust as compared to more evolved stage 
where it has cleared all surrounding matter. A rarer 
environment in a more evolved stage could lead to an 
efficient acceleration of winds via line driving to 1000 
kms -1 . Hence, for more massive and hotter young stars, 
impact due stellar winds could play a vital role on the 
dynamical evolution of disk winds. 

6.4.3. Lack of knowledge of system parameters 

The inner regions around massive young stars are not 
accesible with present day telescopes. Thus, many phys- 
ical parameters for these inner regions are not observa- 
tionally constraint (for e.g the inner disk radius, mass 
density, magnetic field strength) In the present model, we 
follow the notion that high-mass stars form as a scaled- 
up low-mass stars. Thus most of the parameters used 
for the present modeling are derived from estimates for 
low-mass stars. 

6.4.4. Ideal MHD and line force parameters 

Simulations presented here are done using ideal MHD 
approximation. This is in principle fine as ionizations 
fractions are usually high enough to provide a good cou- 
pling between matter and field. However, for simplicity 
we have assumed a constant ionization fraction through- 
out the wind. The radiative force parameters do have an 
explicit dependence on the ionization fraction. The vary- 
ing temperature across the jet implies to varying degree 
of ionization which would modify the line force parame- 
ters and thus the value of the force. In this model, we use 
an idealized assumption of a constant ionization fraction 
throughout and fix the radiative force parameters for all 
simulation runs. While in case of non-magnetic CVs, 
Pereyra & Kallman (2003) has shown that the qualita- 
tive effect of the line force is not altered by incorporat- 
ing local ionization effects for the line force parameters. 
However, the importance of such local ionization effects 
for the case of young massive stars is not known in the 
literature. 

6.4.5. Prescription of the disk dynamics 

A self-consistent modeling would need to incorporate 
the time-evolution of the disk structure in the simulation, 
and to treat accretion and ejection processes simultane- 
ously. Such simulations are currently performed in gen- 
eral applications of jet launching (e.g Casse & Keppens 
2002, Zanni et al. 2007, Murphy et al. 2010), however, it 
is too early to be applicable to massive young stars - one 
reason being the lack of knowledge of the accretion disk 
parameters (e.g. the question whether there is a thin or 
thick disk), another one that also radiative effects would 
have to be implemented in these simulations. 

7. CONCLUSIONS 

We have studied the impact of a line-driven radia- 
tion forces on the acceleration and collimation of a MHD 
jets, around young massive stellar object. Our main mo- 
tivation was to give a solid theoretical understanding for 
the outflow evolution hypothesis presented by Beuther & 
Shepherd (2005). For the radiation forces we have con- 
sidered stellar and disk luminosity. Our basic approach 
was to initially launch a MHD jet from the underlying 
disk surface, and then, after the pure MHD outflow has 



achieved a steady state, to switch on the radiative forces 
and study their effect on the existing MHD outflow. 

We performed a number of simulations with the line 
driving force exerted by stellar radiation. These simu- 
lations were performed for a wide range of physical pa- 
rameters, as i) the central stellar mass, ii) the magnetic 
flux, and iii) the line- force parameter a. In order to 
apply our method of calculating the line driving force 
(CAK approach), we have modified the numerical code 
PLUTO to incorporate appropriate projections of gradi- 
ents of the 2 D velocity field along the light path. Addi- 
tionally, we have consistently implemented these projec- 
tions for different radiation sources properly considering 
the geometry of the radiation field. All these simulations 
have 'floating' and casually consistent inflow boundary 
in which the mass flux is not fixed and is derived by the 
numerical solution. 

Our main conclusions from this analysis are as follows. 

The line driven force from the central star for the pa- 
rameters considered does play a significant role in mod- 
ifying the dynamics in terms of collimation and acceler- 
ation of the outflow. We find that the outflow velocity 
is increased by a factor of 1.5 - 2 by radiation forces as 
compared to the pure MHD flow. Also, the degree of col- 
limation is lowered, visible e.g. in a 30%-change in the 
magnetic flux profile, or the wider opening angle of the 
magnetic field lines. 

Investigating different stellar masses we determine the 
amount of de-collimation by measuring the opening an- 
gles of a typical field lines (that with the highest mass 
flux) at the Aflven and fast magneto-sonic point. We find 
that for a stellar mass of 2OM , the opening angles are 
20° and 16°, respectively. For a 6OM star these values 
increase to 32° and 29°, indicating substantial amount 
of de-collimation due to the increased stellar luminosity. 
This findings confirm the observed evolutionary picture 
for massive outflows in which more massive stars tend 
to have outflows of lower collimation degree. Note that 
for massive young stars, our results do not only consti- 
tute a relation of different stellar masses, but correspond 
also to a time-evolution of outflow systems, as the central 
mass can be substantially increased during massive star 
formation. 

We have also performed simulations with injection of 
fixed mass flux from the disk boundary for various stellar 
masses. We find that the ratio of vertical to radial mass 
flux in these runs decreases from 0.92 to 0.80 with in- 
crease in stellar mass from 25 M to 50 M . This clearly 
indicates the fact that the line driving force from cen- 
tral star plays an influential role in the physical under- 
standing of observed evolutionary picture pertaining to 
outflows from young massive stars. 

So far, the magnetic field strength in the jet forma- 
tion region close-by a massive young star forming is 
an unknown quantity. We therefore have carried out 
simulations with different field strength (plasma- /?o ~ 
1.0,3.0,5.0). We find that for the flow with high mag- 
netic flux the radiative forces do not really affect the colli- 
mation degree. However, for the flow with low magnetic 
flux, the dynamical effect of radiative forces is greatly 
increased. 

We further find that the line force parameter a is very 
critical in determining the magnitude of the line-driven 
forces. Lower values of a leads to an efficient radiative 
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force from the central star and thus decollimate the flow 
to a larger extent as compared to higher a values. Since 
the radiation forces also affect the mass outflow rates for 
simulations, even small change in a may lead to signifi- 
cant changes in mass flux up to ~ 28%. 

Line forces due to the hot accretion disk do not play a 
significant role in controlling the dynamics of the MHD 
outflow, simply because they are orders of magnitude 
smaller then all other forces that affect the flow dynam- 
ically. Implementing high disk radiation forces has (by 
taking into account the inner hotter part of the disk), 
however, shown us that the disk radiation forces will af- 
fect primarily the flow acceleration, and not the flow col- 
limation. 

The source terms for the line-driven forces from the 
star and the disk implemented in our simulations de- 
pends on certain scaling parameters. We find that the 
physical scaling of the jet density is a leading parameter 
that affects the flow dynamics. Large densities make it 
difficult for the line-driving to act efficiently, resulting in 
a flow which is mostly dominated by MHD forces. How- 
ever, less dense inner regions around massive young stars 
would allow efficient radiative line-driving, and thus ac- 
celerate and de-collimate the flow with great effect. 

This paper has addressed a complex problem of jet 
launching from young massive stars. In doing so, we 
have applied a simplified prescription of the radiative 
force rather than a complete radiative transfer calcula- 
tion. One limitation of our model is the lack of observa- 



tional knowledge about various parameters particularly 
very close to the central star. In addition, we have not 
included effects from stellar winds which might exist dur- 
ing the later stages of young massive star evolution. In 
spite of these limitations, we find clear evidence of accel- 
eration and de-collimation of jets launched from massive 
young stars. 

In summary, among all the radiative sources consid- 
ered to study the dynamics of outflow launched from the 
young massive star, we see that the line force from the 
central luminous star is very efficient in de-collimating 
and accelerating the flow. The line force from the under- 
lying disk is not as significant as compared to the stellar 
force. Also, dynamical effects on the outflow due to the 
optically thin electron scattering continuum force from 
the central star and the disk is not significant. Further- 
more, we confirm the observed trend of de-collimation 
seen in outflows from massive stars at different evolu- 
tionary stages. 
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APPENDIX 

BASICS OF CASTOR, ABBOTT AND KLEIN THEORY 
According to the Castor, Abbott and Klein theory (Castor et al. 1975), the force multiplier can be expressed as, 
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where Az^d the Doppler shift, F v the radiation flux at frequency z/, and F the total integrated flux. The optical depth 
parameter T is related to the gradient of velocity, the density in the wind and the ion thermal velocity v t h, 
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The optical depth parameter T can be related to the optical depth of a particular line tl (n) = r]T. The line strength 
r] is the ratio of the line opacity to the electron scattering opacity cr e , while n is the unit vector along the line of 
sight (l.o.s.). The classification of lines in optically thin and thick lines is done on the basis of interaction probabilities. 
Optically thick lines are those which have interaction probability of unity. Lines with optical depths tl < 1 are 
optically thin, and their probability of interaction is tl. Based on this approximation, the force multiplier is separable 
and can be written differently for optical depths very high or very low. When the gas is optically thick, tl > 1, the 
force multiplier depends only on the local dynamical quantities of the flow, 
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while for the optically thin case tl < 1 the force multiplier is independent of the local dynamics, but depends on the 
line strength of individual thin lines, 

M thin (T)= £ f^V). (A4) 
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In general, for a gas distribution with a mixture of optically thick and thin lines, the empirical form of the total 
force multiplier integrated over all lines can be expressed as a power law, M(T) ~ kT~ a , where k and a are line force 
parameters. 
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LINE DRIVING FORCE DUE TO DISK ALONE 

In addition to the line radiative forces from the central star, we also take into account line forces from the underlying 
disk. The underlying disk cannot be considered as a point source but rather is an extended cylindrical source of 
radiation. Further, the disk luminosity varies with the radial distance from the central star. In the present simulations, 
we consider that the underlying disk as a steady-state standard thin disk with a temperature profile given by Shakura 
& Sunyaev (1973). Thus, the energy radiated per unit area D(r) at cylindrical radius r is 
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where M acc is the steady-state accretion rate onto a central star of mass M*. The inner launching radius is lo. In case 
of line forces from the accretion disk, we calculate the radial and vertical radiation flux from the standard disk, S r and 
5 Z , similar to Pereyra et al. (2000). Further, the l.o.s. velocity gradient in the force multiplier applied in case of disks 
is for simplicity reduced to dv z /dz, since the bulk of the radiation flux from the disk is in vertical direction, 
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Here, S r and S z are the radial and vertical flux components. Both depend on the disk luminosity and are implemented 
in the code in their dimensionless form. 

In order to reduce to a dimensionless form, all the physical quantities should be written as a product of its value in 
code units with appropriate scale factor (see Eq. 12). S r and S z are radial and vertical components of the flux emitted 
from the disk surface. (Hachiya et al. 1998, Pereyra et al. 2000) - 
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where D(r 7 ) is the amount of energy radiated per unit area from the standard thin disk (Eq Bl). In the present 
simulations, the accretion disk is treated as a boundary and the accretion of matter in the disk is not considered. The 
mass accretion rate that appears in the above flux radiative components of the disk has to be prescribed on basis of 
dimensionless parameters. 

Macc = ^r e A M . (B5) 



The radiation flux given above has to be written in dimensionless form to be incorporated in the simulations. These 
flux components in the code units are given as follows - 
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Using the above formulations, the dimensionless radial and vertical component of the line force from the disk can 
then be obtained. Their respective contours are shown in Fig. B. 



/^tk^c^c) = /wdiskAGM./O = M c (T)S T , code , 
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Figure Bl. Contours of the line force from disk radiation in the radial (left panel) and vertical (right panel) direction. The contour 
levels are given in physical units. The parameters used are : Qo = 1400.0, a = 0.55, M* = 3OM0, Iq = 0.1 AU, po = 5.0 x 10 -14 gem -3 , 
T e = 0.2369, A = 0.4969, [i = 0.4644,/3 = 5.0 

The force multiplier for the case of disk force in code units is similar to that used for stellar line force (Eq 13). 
However, since we assume that bulk of the photons from the disk move along the vertical z axis, the line of sight 
velocity gradient is approximated to be just due to vertical velocity, 
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The contours of the line force due to the underlying disk in the radial and the vertical direction are shown in figure 
below. 
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